The Stability of an Expanding Circular Cavity and the Failure of Amorphous Solids 
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Recently, the existence and properties of unbounded cavity modes, resulting in extensive plastic 
deformation failure of two-dimensional sheets of amorphous media, were discussed in the context 
of the athermal Shear- Transformation-Zones (STZ) theory. These modes pertain to perfect circular 
symmetry of the cavity and the stress conditions. In this paper we study the shape stability of the 
expanding circular cavity against perturbations, in both the unbounded and the bounded growth 
regimes (for the latter the unperturbed theory predicts no catastrophic failure). Since the unper- 
turbed reference state is time dependent, the linear stability theory cannot be cast into standard 
time-independent eigenvalue analysis. The main results of our study are: (i) sufficiently small per- 
turbations are stable, (ii) larger perturbations within the formal linear decomposition may lead to an 
instability; this dependence on the magnitude of the perturbations in the linear analysis is a result of 
the non-stationarity of the growth, (iii) the stability of the circular cavity is particularly sensitive to 
perturbations in the effective disorder temperature; in this context we highlight the role of the rate 
sensitivity of the limiting value of this effective temperature. Finally we point to the consequences 
of the form of the stress-dependence of the rate of STZ transitions. The present analysis indicates 
the importance of nonlinear effects that were not taken into account yet. Furthermore, the analysis 
suggests that details of the constitutive relations appearing in the theory can be constrained by the 
modes of macroscopic failure in these amorphous systems. 
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I. INTRODUCTION 

Some of the theoretically most fascinating aspects of 
crack propagation in amorphous materials are the insta- 
bilities that are observed in well controlled laboratory 
experiments [l|. Besides some exceptions, (see for exam- 
ple [2, S \M ^^^ s-lso 01) J it would be fair to say that 
the observed instabilities are still poorly understood. It 
is the opinion of the present authors that the reason for 
the relative lack of understanding is that the theory of 
crack propagation did not treat cracks as moving free 
boundaries whose instabilities stem from the dynamics 
of the free boundary itself. Instead, "crack tip dynam- 
ics" were replaced by energy balance within the theory 
of Linear Elastic Fracture Mechanics Q, together with 
an ad-hoc "law" of one nature or another as to where a 
crack is supposed to move. 

In principle this undesirable state of affairs can be 
greatly improved within the Shear- Transformation-Zones 
(STZ) theory of amorphous materials [1, i, E, [HI ■ This 
theory treats developing cracks or growing cavities as 
free boundaries of a material in which both elasticity 
and plasticity are taken into account, preserving all the 
symmetries and conservation laws that promise a pos- 
sibly correct theory of amorphous materials driven out 
of mechanical equilibrium. This theory in its various ap- 
pearances was compared to a number of experiments and 
simulations (sec below), with a growing confidence that 
although not final, STZ theory is developing in the right 
direction. Indeed, the application of a highly simplified 
version of STZ theory to crack propagation resulted in 
physically interesting predictions, explaining how plas- 
ticity can intervene in blunting a crack tip and result- 
ing in velocity selection [l2| . The application of the full 



fledged theory of STZ to crack propagation is still daunt- 
ing (although not impossible) due to the tensorial nature 
of the theory and the need to deal with an extremely stiff 
set of partial differential equations with a wide range of 
time-scales and length-scales involved. For that reason 
it seemed advantageous to apply the full theory to a sit- 
uation in which the symmetries reduce the problem to 
an effectively scalar theory; this is the problem of a cir- 
cular cavity developing under circular symmetric stress 
boundary conditions [ij, llJI ■ While this problem docs 
not reach the extreme conditions of stress concentration 
that characterizes a running slender crack, it still raises 
many physical issues that appear also in cracks, in par- 
ticular the give-and-take between elasticity and plastic- 
ity, the way stresses are transmitted to moving bound- 
aries (in apparent excess of the material yield stress) and 
most importantly for this paper, the possible existence 
of dynamical instabilities of the moving free boundary. 
This last issue might also be connected to the difference 
between ductile and brittle behaviors. In the former, a 
growing cavity is likely to remain rather smooth, whereas 
in the latter, one may expect an instability resulting in 
the growth of "fingers" , possibly ending up being cracks. 
It is one of the challenges of the present paper to ex- 
amine whether the theory may predict a transition, as a 
function of material parameters or a constitutive relation, 
between these two types of behavior. 

Note that we have chosen to study the problem 
in a purely 2-dimensional geometry; recently quasi 2- 
dimensional systems exhibited interesting failure dynam- 
ics in laboratory experiments, where the 3'rd dimension 
appears irrelevant for the observed phenomena [5J, |l5| . 
Our motivation here is however theoretical, to reduce the 
unnecessary analytic and numerical complications to a 



iTLinimum and to gain insight as to the main physical ef- 
fects under the assumption that the thin third dimension 
in real systems does not induce a catastrophic change in 
behavior. When this assumption fails, as it does in some 
examples c.f. |l|, the analysis must be extended to in- 
clude the third dimension. This is beyond the scope of 
this paper. 

In Sec. |TT]we present the equations that describe the 
problem at hand and specify their boundary conditions. 
Particular attention is paid to distinguish between the 
general Eulerian formulation which is model-independent 
(Subsec. Ill A[) and the constitutive relations involving 
plasticity where the STZ model is explained (Subsec. 
IIIBj) . This section finishes with the presentation of the 
unperturbed problem, preparing the stage for the linear 
stability analysis which is discussed in Sec. IIIII In this 
section we present a general analysis where inertia and 
elastic compressibility effects are taken into account. In 
Appendix |B] we complement the analysis by consider- 
ing the "quasistatic" (when the velocity of the bound- 
ary is sufficiently small) and incompressible case (when 
the bulk modulus is sufficiently large) and show that 
both formulations agree with one another in the relevant 
range. The results of the stability analysis are described 
in detail in Sec. IIVI and a few concluding remarks are 
offered in Sec. IVl 



II. EQUATIONS AND BOUNDARY 
CONDITIONS 



A. General formulation 

We start by writing down the full set of equations for 
a general two-dimensional elasto-viscoplastic material. A 
basic assertion of this theory is that plastic strain tensors 
in such materials are not state variables since their values 
depend on the entire history of deformation. Thus, one 
begins by introducing the total rate of deformation tensor 
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where v{r,t) is the material velocity at the location r 
at time t and T denotes here the transpose of a tensor. 
This type of Eulerian formulation has the enormous ad- 
vantage that it disposes of any reference state, allowing 
free discussion of small or large deformations. As is re- 
quired in an Eulerian frame we employ the full material 
time derivative for a tensor T, 
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'Vt ^ ~dt 
where u) is the spin tensor 
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For a scalar or vector quantity V the commutation 
with the spin tensor vanishes identically. The Eulerian 
approach allows a natural formulation of moving free 



boundary problems; this will be shown to lead to a signifi- 
cant advance compared to more conventional treatments. 
The plastic rate of deformation tensor I?p' is intro- 
duced by assuming that the total rate of deformation 
tensor i)'°* can be written as a sum of a linear elastic 
and plastic contributions 
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We further assume that D^' is a traceless tensor, corre- 
sponding to incompressible plasticity. All possible ma- 
terial compressibility effects in our theory are carried by 
the elastic component of the deformation. The compo- 
nents of the linear elastic strain tensor e'^' are related to 
the components of stress tensor, whose general form is 
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where K and /i are the two dimensional bulk and shear 
moduli respectively. The tensor s is referred to hereafter 
as the "deviatoric stress tensor" and p as the pressure. 
The equations of motion for the velocity and density are 
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In order to prepare the general set of equations for the 
analysis of a circular cavity we rewrite the equations in 
polar coordinates. For that aim we write 
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where e^ and eg arc unit vectors in the radial and az- 
imuthal directions respectively. These expressions en- 
able us to represent the divergence operator V- in the 
equations of motion and the covariant derivative v ■ V 
in the material time derivative of vectors and tensors. 
Some care should be taken in evaluating these differen- 
tial operators in polar coordinates since the unit vectors 
themselves vary under differentiation according to 

drBr = 0, drBe = 0, deSr = Bg, dgee = -Br . 

(2.10) 
We then denote Srr = —s, sgg = s, Sro — sgr = T and 
using Eqs. (|2.5|) we obtain 
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In this notation the equations of motion (|2.7p can be 
rewritten explicitly as 
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where V • <t is calculated explicitly in Appendix [X] 

Equations (j2.4p can be rewritten in components form 
as 
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Here the components of the total rate of deformation 
tensor are related to the velocity according to Eqs. (|2.ip 
as 
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where the components of the spin tensor oj in Eq. 
are given by 
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The calculation of the tensor tJ-Ve*^' is presented in Ap- 
pendix |A1 the linear elastic strain components of Eqs. 
p.6p are given by 
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Since most of the materials of interest have a large 
bulk modulus K, i.e. they are almost incompressible, we 
assume that the density is constant in space and time 
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Therefore, Eq. (|2.8p is omitted. Finally, the existence of 
a free boundary is introduced as the following boundary 
conditions 

a.jHj = , (2.18) 

where n is the unit normal vector at the free boundary. 



B. Viscoplastic constitutive equations: 
The athermal STZ theory 



Up to now we have considered mainly symmetries and 
conservation laws. A general theoretical framework for 
the elasto-viscoplastic deformation dynamics of amor- 
phous solids should be supplemented with constitutive 
equations relating the plastic rate of deformation tensor 
dp'' to the stress and possibly to other internal state 
fields. We use the constitutive equations of the re- 
cently proposed athermal Shear Transformation Zones 
(STZ) theory [ll|. This theory is based on identify- 
ing the internal state fields that control plastic deforma- 
tion. The basic observation is that stressing a disordered 
solid results in localized reorganizations of groups of par- 
ticles. These reorganizations occur upon surpassing a 
local shear threshold, and when they involve a finite ir- 
reversible shear in a given direction, we refer to them as 
an "STZ transition". Once transformed, due to a local 
redistribution of stresses, the same local region resists 
further deformation in that direction, but is particulary 
sensitive to shearing transformation if the local applied 
stress reverses its direction. Thus an STZ transition is 
conceived as a deformation unit that can undergo config- 
urational rearrangements in response to driving forces. 
Furthermore, the stress redistribution that accompanies 
an STZ transition can induce the creation and annihi- 
lation of other local particle arrangements that can un- 
dergo further localized transitions; these arrangements 
are formed or annihilated at a rate proportional to the 
local energy dissipation (recall that thermal fiuctuations 
arc assumed to be absent or negligible). In this sense 
the interesting localized events need not depend on "pre- 
existing" defects in the material, but can appear and dis- 
appear dynamically in a manner that we describe math- 
ematically next. 

This picture is cast into a mathematical form in terms 
of a scalar field A that represents the normalized density 
of regions that can undergo STZ transitions, a tensor m 
that represents the difference between the density of re- 
gions that can undergo a transition under a given stress 
and the reversed one, and an effective disorder temper- 
ature X that characterizes the state of configurational 
disorder of the solid [lg| . The present state of the theory 
relates these internal state fields, along with the devia- 
toric stress tensor s, to the plastic rate of deformation 



tensor D^' according to 
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This equation represents the dependence of the plastic 
rate of deformation on the current stress Sij and the re- 
cent history encoded by the internal state tensorial field 
m. This field acts as a back-stress, effectively reduc- 
ing the local driving force for STZ transitions, up to the 
possible state of jamming when the whole parentheses 
vanishes. The parentheses provides information about 
the orientation of the plastic deformation. The function 
C{s) determines the magnitude of the effect, and is re- 
discussed below. The field A appears multiplicatively 
since the rate of plastic deformation must be proportional 
to the density of STZ. The second equation describes the 
dynamics of the internal back stress field 
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This equation captures the dynamical exchange of stabil- 
ity when the material yields to the applied stress. The 
equation has a jammed fixed point when the plastic de- 
formation vanishes, in agreement with the state of STZ 
being all in one orientation, without the production of a 
sufficient number of new ones in the other orientation. 
The jammed state is realized when the applied stress 
is below the yield stress. When the stress exceeds the 
threshold value the stable fixed point of this equation 
corresponds to a solution with non- vanishing plastic rate 
of deformation. This state corresponds to a situation 
where enough STZ are being created per unit time to 
allow a persistent plastic fiow. The quantity F repre- 
sents the rate of STZ production in response to the fiow 
D^K The next equation, for the density of STZ A, is an 
elementary fixed point equation reading 
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The unique fixed point of this equation is the equilibrium 
solution A = e~^l^ where x is a normalized temperature- 
like field which is not necessarily the bath temperature 
when the system is out of thermal and/or mechanical 
equilibrium. The last equation is for this variable, read- 
ing 

Toco^ = eoAr(sy,my) [xoo {tqD'^^) - x] , 



with DP' = 
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This is a heat-like equation for the configurational de- 
grees of freedom; it is discussed in detail below. Here and 
elsewhere we assume that quantities of stress dimension 



are always normalized by the yield stress Sy\ this is justi- 
fied as the STZ equations exhibit an exchange of dynamic 
stability from jamming to flow at 5 = 1, i.e. at a stress 
that equals to Sy ^. The set of Eqs. ^I^-^I^ is a 
tensorial generalization of the effectively scalar equations 
derived in [lll | ; such a generalization can be obtained by 
following the procedure described in Ref. ^M,- In these 
equations, tq is the elementary time scale of plasticity, 
eg is a dimensionless constant and cq is a specific heat in 
units of kB per particle. 

A weak point of the theory is the lack of a first- 
principle derivation that determines the function C{s) in 
Eq. (|2.19p . which lumps together much of the micro- 
scopic physics that controls the stress-dependent rate of 
STZ transitions. Our theory constrains it to be a sym- 
metric function of s that vanishes with vanishing deriva- 
tives at s = 0, due to the athermal condition that states 
that no transitions can occur in a direction opposite to 
the direction of s [ll|]. This constraint is not sufficient, 
however, to determine C[s). To appreciate the uncertain- 
ties, recall that STZ transitions are relaxation events, 
where energy and stress arc expected to re-distribute. 
Even without external mechanical forcing, aging in glassy 
systems involves relaxation events that are poorly under- 
stood [l3|- The situation is even more uncertain when 
we deal with dynamics far from mechanical equilibrium. 
The best one can do at present is to choose the function 
C{s) by examining its infiuence on the resulting macro- 
scopic behaviors [l9| . Thus in this paper we will examine 
the sensitivity of the stability of the expanding cavity to 
two different choices of C{s). At present we use the one- 
parameter family of functions, C{s) = J-{s\Q^ proposed 
in [nl 



n^\0^ 



^CH 



C! 



s| 



|s|-SQ)si exp(-Csa)(isa . (2.23) 



The integral is over a distribution of transition thresholds 
whose width is controlled by a parameter C (and see [11 1 
for details). For finite values of C there can be nonzero 
sub-yield plastic deformation for | s | < 1 . This behavior is 
well documented in the literature cf. [20| in the context 
of experimental stress-strain relations and plastic defor- 
mations. We note that for s very small or very large, 
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In Sec. lIVCI wc propose a different one-parameter family 
of functions Q{s] A) and study in detail the implications 
of this different choice on the stability of the expanding 
cavity. 

Eq. (|2.22p deserves special attention. It is a heat-like 
equation for the effective disorder temperature x with a 
fixed-point Xoo which is attained under steady state de- 
formation. This refiects the observations of Ref. [2l| . 
where the effective temperature x "^as shown to attain 
a unique value in the limit t^D^^ — > 0, where to was the 
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FIG. 1: (Color online) A typical relation between Xoo and 
logio(ro-D'''j for a temperature significantly smaller than the 
glass transition temperature. Data courtesy of T. Haxton and 
A. Liu. Note that the data were scaled properly. 



particles vibrational time scale. Indeed, in most applica- 
tions, realistically imposed inverse strain rates are much 
larger than the elementary time scale to, i.e. io-D^' ^ 1. If 
we identify our tq with the vibrational time scale io (see 
for example |22l | ) , we conclude that Xcxi can be taken as a 
constant, independent of the plastic rate of deformation. 
This assumption was adopted in all previous versions of 
STZ theory. Note also that a low plastic rate of deforma- 
tion is associated with s—>l+, i.e. a deviatoric stress that 
approaches the yield stress from above. However, the 
situation might be very different in free boundary evolu- 
tion problems, where high stresses concentrate near the 
boundary, reaching levels of a few times the yield stress. 
Estimating x in the typical range of 0.1 — 0.15 [23, l23| . 
e~^l^ is in the range 10"'* — 10"'^. Therefore, estimating 
the other factors in Eq. (|2.19p . for the high stresses near 
the free boundary, in the range 1 — 10, we conclude that 
tqD'p^ can reach values in the range 10^^—10^^. Very re- 
cent simulations [2J] demonstrated convincingly that in 
this range of normalized plastic rates of deformation, Xao 
shows a considerable dependence on this rate, see Fig. [T] 
Since x affects plastic deformation through an exponen- 
tial Boltzmann-like factor, even small changes of Xoa in 
Eq. (|2.22p can generate significant effects [2^ . This issue 
is of particular importance for the question of stability 
(and localization) under study since the strain rate sensi- 
tivity of Xoo niight incorporate an instability mechanism; 
fluctuations in the plastic rate of deformation, caused for 
example by fluctuations in Xi can induce, through Xoo, 
a further localized increase in plastic deformation and so 
on. This intuitive idea will be studied in the analysis to 
follow. 

The set of Eqs. ((2J9l) - ((2:221) 01 (and slight variants) 
was shown to capture viscoelastic behavior in a variety 
of examples. These include small stress and finite plas- 
ticity at intermediate stresses |27| , a transition to flow at 



the yield stress (as discussed above) [ll|, the deforma- 
tion dynamics of simulated amorphous silicon |22l |. the 
necking instability |28|. the deformation dynamics near 
stress concentrations (1311 . the cavitation instability [1J| 
and strain localization [29| . In this work we focus on the 
implications of these constitutive equations on the sta- 
bility of propagating free boundaries in relation to the 
failure modes of amorphous solids. 



C. The unperturbed problem 

In this subsection we adapt the general theory to the 
circular symmetry of the unperturbed expanding cavity 
problem. We consider an infinite medium with a circular 
cavity of radius i?^*'^(i), loaded by a radially symmet- 
ric stress CT°° at infinity. The superscript (0) in all the 
quantities denotes the fact that they correspond to the 
perfectly symmetric case that is going to be perturbed 
later on. For the perfect circular symmetry the velocity 
field v''''(r,6') is purely radial and independent of the 
azimuthal angle d, i.e. 
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This symmetry also implies that 
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(2.26) 
and all the diagonal components are independent of Q. 
Eqs. (|2.4p . after a simple manipulation, can be rewritten 
as 
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where we have defined 
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The equations of motion (|2.12p reduce to 
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The boundary conditions are given by 
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The initial conditions are chosen to agree with the solu- 
tion of the static linear-elastic problem, i.e. 
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This choice reflects the separation of time scales between 
elastic and plastic responses. This separation of time 
scales can be written explicitly in terms of the typical 
elastic wave speed, the radius of the cavity and the time 
scale of plasticity: 
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Finally, the rate of the cavity growth is simply deter- 
mined by 
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For the circular symmetry, the STZ equations (|2.19p - 
([2:20)) reduce to 
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Note that the x a-nd D^' equations contain a factor 
of the small STZ density eoA^^^ , which implies they are 
much stiffer than the m and A equations. Therefore, 
whenever the advection terms can be neglected this sep- 
aration of time scales [13 allows us to replace the equa- 
tions for m^^' and A^°-' by their stationary solutions 
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Note that Eq. (|2.36p has two stable fixed-point solu- 
tions given by Eq. ()2.39p . where we used Eq. ()2.40|) 
and omitted the advection term. The transition between 
these two solutions corresponds to a transition between 
a jammed and a plastically flowing state for a deviatoric 
stress below and above the yield stress respectively [111 . 
Eq. (|2.35p exhibits the corresponding solutions in terms 
of the plastic rate of deformation, zero and finite, below 
and above the yield stress respectively. 

The unperturbed problem was studied in detail in Ref. 
[14| . It was shown that for stresses a°° smaller than a 



threshold value ct*^ ~ 5 the cavity exhibits transient dy- 
namics in which its radius approaches a finite value in a 
finite time. When this happens the material is jammed. 
On the other hand, for a°° > cr^^ the cavity grows with- 
out bound, leading to a catastrophic failure of the ma- 
terial, accompanied by large scale plastic deformations. 
We stress that to our knowledge this mode of failure by 
propagating a plastic solution is new, apparently not re- 
lated to other recently discovered failure fronts [301 ■ O'^s 
major goal of the present study is to analyze the sta- 
bility of the unbounded growth modes that result from 
this cavitation. However, we are also interested in the 
range ct°° < tr*'' where the unperturbed theory predicts 
no catastrophic failure. In this range, a failure can still 
occur if the cavity, prior to jamming, loses its perfect cir- 
cular symmetry in favor of relatively slender propagating 
"fingers" . In that case, stress localization near the tips 
of the propagating "fingers" can lead to failure via frac- 
ture. Such a scenario is typical of brittle fracture where 
the stress localization due to the geometry of the defect 
drives crack propagation that might lead to macroscopic 
failure. 



III. LINEAR STABILITY ANALYSIS 

We derive here a set of equations for the linear per- 
turbations of the perfect circular symmetry where both 
inertia and clastic compressibility effects are taken into 
account. In Appendix [B] we complement the analysis 
by considering the quasi-static and incompressible case. 
This case is mathematically more involved as it contains 
no explicit time evolution equation for the velocity and 
the pressure fields. By comparing the results of the two 
formulations we test for consistency and obtain some de- 
gree of confidence in the derivation and the numerical 
implementation of the equations presented in this sec- 
tion. 



A. Equations of motion and kinematics 

The quantities involved in the problem are the tensors 



— S T 
T S 



DPi 



Dl', DP' 



(3.1) 



as well as the pressure p{r,t), the velocity v{r,t) and 
the location of the free boundary R{9,t). We start by 
expanding all these quantities as follows 

R{e,t) = i?('')(t) + e"^''i?(i'(i) , 

s{r,0,t) = .s(°)(r,i)+e*"^(i)(r,t) , 

T{r,e,t) = ie"'''T^^^r,t) 

p{r,e,t) = p^°'>{r,t) + e"'''p^^'>{r,t) , 

vg{r,0,t) = ie"^%i'\r,t) , 



;(0) 



(1), 



DP\r,e,t) = D'P^^\r,t) + e'""^ DP^^\r,t) 



D%{r,0,t) 



AnO 



D 



pi (1) 



{r.t) 



(3.2) 



Here all the quantities with the superscript (1) are as- 
sumed to be much smaller than their (0) counterparts 
and n is the discrete azimuthal wave-number of the per- 
turbations. The small perturbation hypothesis results in 
a formal linear decomposition in which each linear mode 
of wave- number n is decoupled from all the other modes. 
When nonlinear contributions arc non-negligible, all the 
modes become coupled and the formal linear decomposi- 
tion is invalid. 

We expand then the equations of motion (|2.12p to first 
order to obtain 
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We proceed by expanding Eqs. (|2.4[) to first order, which 
after a simple manipulation yields 
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At this point wc derive an evolution equation for 
the dimensionlcss amplitude of the shape perturbation 
rW/r("). To that aim we note that 



R = Vr{R)+0 






(3.8) 



Expanding this relation using Eqs. p.2p . we obtain to 
zeroth order Eq. (|2.34p and to first order 



i?(i)(t) 



,(l)rR(") 



(i?("') + i?« 



dvrjR^"^) 

dr 



(3.9) 



Therefore, we obtain 
dt \R(o)J ^ 
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This is an important equation since a linear instability 
manifests itself as a significant increase in i?(^)/i?'^°) such 
that nonlinear terms become non-negligible. Note that 
the two last terms in the square brackets are always neg- 
ative, therefore an instability can occur only if the first 
term in the square brackets is positive with absolute value 
larger than the sum of the two negative terms. Moreover, 
recall that the problem is non-stationary, implying that 
all the zeroth order quantities depend on time. 

In order to derive the boundary conditions for the com- 
ponents of the stress tensor field we expand to linear or- 
der the normal unit vector n (not to be confused with 
the discrete wave-number n) and tangential unit vector 
t at the free boundary, obtaining 



n = 1, — i 



.i?(i 



-ne 



inO 



i?(i 



?7,e™^l 



Eqs. (|2.18p . expanded to first order, translate to 
s(i)(i?("')+pW(i?(°)) = 



9s(o)(i?(")) ap(o)(i?(o)) 



dr 



dr 



l(i?("')=7lLs(")(i?("))-p(")(i?("))' 



i?(l) 
i?(")' 



(3.11) 

(3.12) 
(3.13) 



In addition, all the first order fields decay as r —>■ oo. 
The initial conditions are determined by the perturbation 
scheme that is being studied. 

To avoid dealing with an infinite and time-dependent 
domain we applied the following time-dependent coordi- 
nate transformation 



e = R{t)/r 



(3.14) 



This transformation allows us to integrate the equations 
in the time- independent finite domain ^G [0, 1], with the 
price of introducing new terms in the equations. Con- 
trolling the equations at small distances required the in- 
troduction of an artificial viscosity on the right-hand-side 
(RHS) of Eq. (HJI). The term introduced is pif7^, with 
7] chosen of the order of the square of space discretiza- 
tion over the time discretization. This introduces zeroth 
order contributions on the RHS of Eq. (|2.30p and first 
order contributions on the RHS of Eqs. ((a!3)) - ((0)) . 



B. 



Linear perturbation analysis of the STZ 
equations 



The only missing piece in our formulation is the per- 
turbation of the tensorial STZ equations. In addition to 



the fields considered up to now, the analysis of the STZ 
equations includes also the internal state fields 
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Therefore, in addition to Eqs. (|3.2p we have 

m{r,e,t) = m(°)(r,<)+e"V(i)(r,t) , 
mre{r,e,t) = le'-^'m^^^ {r,t) , 



X(r, e, t) = x'"' (r, t) + e™V^'^ (r, i) • (3.16) 



(3.15) Accordingly we expand C{s) (assuming s'"' > 0) in the 
form 



C(s) = C(sW + e™^(i)) ~ C(sW) + ^ (5(0)) e-''s(i), 

(3.18) 
where 



We then expand systematically Eqs. (|2.19|) - (|2.23|1 . 
First, we have 
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Substituting the last three equations into (|2.19p and ex- 
(3.17) panding to first order, we obtain 
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We then expand F in the form 

r = r(") + e^'^Fd) with F(") = ^^'' ,-^, , 
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Eq. p.20p is now used to obtain 
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Using Eq. (|2.2ip we obtain 
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Expanding 1?^', similarly to Eq. (|3.17p . we obtain 
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Accordingly we expand Xoo (tq^'^') (with 13^' > 0) in 
the form 

Xoo (ro^^') = Xoo (roi^P'^"^ + e^"^oi5^'^'^) (3.27) 



Then, using Eq. (|2.22p we obtain 
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Thus, Eqs. (13:201) - (13:21]) ■ (13:231) - (13:241) . (13:251) and (13:281) 
constitute our equations for the dynamics of the first or- 
der STZ quantities. 

These equations already reveal some interesting fea- 
tures. First note that the coupling between D^' (which 
is the quantity that is expected to be of major importance 
in determining Vr in Eq. (|3.10p through Eqs. (|3.5p - 
1^) and x^^\ m'^^^ depends on C{s'^°^). This means 
that the strength of the coupling depends on (. Sim- 
ilarly, the coupling between D^^ and s'^^ depends on 
dC (s^"-*) /ds which is also a function of (. These observa- 
tions demonstrate the importance of the precise form of 
the function C{s). This issue is further discussed in Sec. 
IIVCI Finally, note that whenever the advection terms 
can be neglected, the known separation of time scales 
[ill allows us to use Eqs. ([1311) -(EHOl) and to replace 
the equations for m^^' , m^g and A'^' by their stationary 
solutions 
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FIG. 2: (Color online) Upper panel: The ratio i?^^'/^'"' as 
a function of time for n — A and a°° =2.1,4.1 and 6.1. Note 
that the last value is above the cavitation threshold. Lower 
panel: The ratio 7?'^^^'"' as a function of time for cr°° =4.1 
and n = 2, 4 and 8. 
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In the next section we summarize the results of our 
analysis of the equations derived in Sec. IlICi IIII Al and 

mm 



IV. RESULTS 

We arc now ready to present and discuss the results 
of the stability analysis of the expanding circular cavity. 
The full set of equations was solved numerically as dis- 
cussed above. Time and length are measured in units of 
To and i?*^°^(i = 0) respectively. A and m are set initially 
to their respective fixed-points. The material-specific pa- 



rameters used are eg 

p = i, x(°)=o.ii,x 



1, CO^l, fl/Sy^bO, K/Sy 



100, 



1; X^°^ =0.11, Xoo = 0.13 and C = 7, unless otherwise 



stated. In Subsec. IIVAI we study perturbations of the 
shape of the cavity and of the effective temperature x- In 
Subsec. IIV Bl we study the effect of the rate dependence 
of Xoo on the stability analysis and in Subsec. IIVCI we 
analyze the effect of the stress-dependent rate function 
ds). 



A. Perturbing the shape and x 

Studying the linear stability of the expanding cavity 
can be done by selecting which fields are perturbed and 
which are left alone. In practice each of the fields involved 
in the problem may experience simultaneous fluctuations, 
including the radius of the cavity itself. Therefore, one 
of our tasks is to determine which of the possible pertur- 
bation leads to a linear instability. To start, we perturb 
the radius of the expanding cavity at f = while all the 
other fields are left alone. In Fig. [2] we show the ratio 
i?^^^/i?'^°-' as a function of time for various loading levels 
a°° (both below and above the cavitation threshold) and 
wave-numbers n. The initial amplitude of the pertur- 
bation was set to i?'^'/i?(°) = 10"'^. The observation is 
that the ratio R^^^ /R^^^^ does not grow in time in all the 
considered cases were the radius was perturbed, implying 
that here the circular cavity is stable against shape per- 
turbations. Note that R^^^ /R^^^^ decays faster for larger 
n and for larger cr°°. Also note that for (t°° = 6.1. i.e. for 
unbounded zeroth order expansion, the ratio i?^^^/i?(°^ 
decays to zero while below the cavitation threshold this 
ratio approaches a finite value. The latter observations 
means that when the material approaches jamming (with 
RS^' attaining a finite value in a finite time) the pertur- 
bation have not yet disappeared entirely. 

We stress at this point the non-stationary nature of the 
problem in which R^^*' (t) is an increasing function of time. 
Thus, even if the absolute magnitude of the amplitude of 
the shape perturbation R^^'{t) increases with time, an 
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FIG. 3: (Color online) i?^^' as a function of time for a°° = 2.1 
and n = 4. 



instability is not automatically implied; i?'^'(t) should 
increase sufficiently faster than R^^\t) in order to imply 
an instability. To exemplify this feature of the problem, 
we present in Fig. [3]i?(i)(t) for a°° = 2.1 and n = 4. It is 
observed that even though R^^^ increases, the smallness 
parameter i?'-'^'/_R^°) decreases, see Fig. [31 Note also 
that i?^^^ does not increase exponentially as expected in 
stationary linear stability analysis, but rather tends to 
asymptote to a constant. 

Next we have tested the stability of the expanding cav- 
ity against initial perturbations in the velocity field or in 
the stress field. The results were quantitatively similar 
to those for the shape perturbations summarized in Figs. 
[H and [31 all implying linear stability. 

In light of these results, we concentrated then on the 
effect of perturbations in the STZ internal state fields. 
Since the dynamics of the tensor m are mainly deter- 
mined by the deviatoric stress field s, we focus on fluctu- 
ations in the effective disorder temperature x- This may 
be the most liable field to cause an instability. Indeed, 
in Ref. [29[ it was shown that x perturbations control 
strain localization in a shear banding instability. Qual- 
itatively, an instability in the form of growing "fingers" 
involves strain localization as well; plastic deformations 
are localized near the leading edges of the propagating 
"fingers". In Ref. 0, based on the data of Ref. [H], 
it was suggested that the typical spatial fluctuations in 
X have an amplitude reaching about 30% of the homo- 
geneous background x- Obviously wc cannot treat such 
large perturbations in a linear analysis and must limit 
ourselves to smaller perturbations. In Fig. [H we show 
the ratio i?'^^'/i?'^*'' as a function of time for a perturba- 
tion of size x'^Vx^°^ = 0.03, introduced at time i = 0. 
The wave- number was set to n = 4 and a°° was set both 
below and above the cavitation threshold. First, note 
that for both loading conditions i?(^'/i?(°) increases on 
a short time scale of about IOOOtq, a qualitatively diffcr- 
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FIG. 4: (Color online) The ratio R^^^/R^°^ as a function of 
time for a perturbation of size x'^Vx''^' =0.03 and n — 4, in- 
troduced at time t = 0. The solid line corresponds to a°° =4.1 
(below the cavitation threshold) and the dashed line corre- 
sponds to cr°° =6.1 (above the cavitation threshold). 



ent behavior compared to the system's response to shape 
perturbations. Second, note the qualitatively different 
response below and above the cavitation threshold. In 
the former case, i?(-'^'/i?'^°^ increases monotonically, ap- 
proaching a constant value when R^^' attains a finite 
value (i.e. jamming). In the latter case, R^^^'/R^^') in- 
creases more rapidly initially, reaches a maximum and 
then decays to in the large t limit. Therefore, in spite 
of the initial growth of R^^yR,^^\ for this magnitude of 
X perturbations, the expanding circular cavity is linearly 
stable; below the cavitation threshold the relative mag- 
nitude of the deviation from a perfect circular symmetry 
i?*^^-'/^*'"'' tends to a finite constant, i.e. a shape per- 
turbation is "locked in" the material, while above the 
threshold the cavity retains its perfect circular symmetry 
in the large t limit. Nevertheless, in light of the signifi- 
cant short time increase in i?(^^/i?(°) (here up to 0.6%), 
we increased the initial (t = 0) x perturbation to the 
range X /x =0.05 — 0.06, in addition to shape per- 
turbations of a typical size of R'-^'^/R'-^'^ = 0.02 - 0.03. 
In these cases i?'^^^/i?'^°) grows above 5%; even more im- 
portantly, the field x'-^H'";^) (^^ ''^ell as other fields in 
the problem) becomes larger than 0.1x^''-'(r, 6*) near the 
boundary of the cavity, invalidating the small perturba- 
tion hypothesis behind the perturbativc expansion and 
signaling a linear instability. Naturally, this breakdown 
of the linearity condition takes place firstly near a peak 
of the ratio i?'^^^/i?(°\ similar to the one observed in Fig. 

SI 

We thus propose that sufficiently large perturbations 
in the shape of the cavity and the effective disorder tem- 
perature Xj but still of formal linear order, may lead to 
an instability. This dependence on the magnitude of the 
perturbations in a linear analysis is a result of the non- 



11 



stationarity of the growth. Another manifestation of the 
non-stationarity is that even in cases where we detect 
an instabihty, it was not of the usual simple exponen- 
tial type where an eigenvalue changes sign as a function 
of some parameter (or group of parameters). Combined 
with the evidence for the existence of large fluctuations 
in X @j [23 1 the present results indicate that it will 
be worthwhile to study the problem by direct boundary 
tracking techniques where the magnitude of the pertur- 
bation is not limited. 

We conclude that the issue of the stability of the ex- 
panding cavity can be subtle. Sufficiently small pertur- 
bations are stable, though there is a qualitative differ- 
ence in the response to perturbations in the effective 
disorder temperature x, where the ratio /J^^^/i?'^"^ in- 
creases (at least temporarily), and other perturbations, 
where R'^^^/R^^'> decays. We have found that for large 
enough x perturbations combined with initial shape per- 
turbations, but still within the formal linear regime, the 
growth of i?(-'^'/i?'^°) takes the system beyond the linear 
regime, making nonlinear effects non-negligible and sig- 
naling an instability. This observation is further sup- 
ported by the existence of large % fluctuations discussed 
in [23, [201 . Note that none of these conclusions depend 
significantly on variations in eo and cq. Moreover, per- 
turbing the expanding cavity at times different than t = 
or introducing a pressure inside the cavity instead of a 
tension at infinity did not change any of the results. 



B. The effect of the rate dependence of Xa 



The analysis of Sec. IIV Al indicates the existence of 
a linear instability as a result of varying the magnitude 
of the perturbations, mainly in %, and not as a result of 
varying material parameters. Here, and in Sec. llVCi 
we aim at studying the effect of material-specific prop- 
erties on the stability of the expanding cavity. Up to 
now we considered Xoo as a constant parameter. How- 
ever, as discussed in detail in Sec. IIIBl the plastic rate 
of deformation near the free boundary can reach val- 
ues in the range where changes in Xoo were observed. 
Therefore, we repeated the calculations using the func- 
tion Xoo(to-D^O plotted in Fig. [T] In Fig. \5\we compare 
R^^^ /R^^^ as a function of time with and without a plas- 
tic rate of deformation dependence of Xoo, both above 
and below the cavitation threshold. The initial pertur- 
bation has x'^^Vx^^'' = 0-03 and n~ 4. Both below and 
above the cavitation threshold the plastic rate of defor- 
mation dependent Xoo{toD^^) induces a stronger growth 
of /j'^^i?'"', though the effect is much more significant 
above the threshold. This is understood as significantly 
higher rate of deformation is developed above the cavita- 
tion threshold, where unbounded growth takes place [ij] , 
compared to below the threshold where the rate of defor- 
mation vanishes at a finite time. We note that the depen- 
dence of Xoo on dp'' affects both the zeroth and first order 
solutions such that _R^°) and i?^^) increase. Our results 
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FIG. 5: (Color online) Upper panel: R^^^/R'-°^ as a func- 
tion of time for a°° = 4.1 (below the cavitation threshold). 
The solid line corresponds to a constant Xoo and the dashed 
line corresponds to the plastic rate of deformation depen- 
dent Xoo(to-D''') of Fig. [1] The initial perturbation has 
X /x ~ 0-03 and n = 4. Lower panel: The same with 
(7°° = 6.1 (above the cavitation threshold). 



show that i?^"'^' is more sensitive to this effect than R^^', 
resulting in a tendency to lose stability at yet smaller 
perturbations. We conclude that the tendency of Xoa to 
increase with the rate of deformation plays an important 
role in the stability of the expanding cavity and might 
be crucial for other strain localization phenomena as the 
shear banding instability [23. Moreover, this material- 
specific dependence of Xodj that was absent in previous 
formulations of STZ theory, might distinguish between 
materials that experience catastrophic failure and those 
that do not, and between materials that fail through a 
cavitation instability [Tj] and those who fail via the prop- 
agation of "fingers" that may evolve into cracks. This 
new aspect of the theory certainly deserves more atten- 
tion in future work. We note in passing that recently an 
alternative equation to Eq. (|2.22p for the time evolution 
of the effective temperature x "^^s proposed in light of 
some available experimental and simulational data [3ll |. 
Preliminary analysis of the new equation in relation to 
the stability analysis performed in this paper indicates 
that the circular cavity does become linearly unstable 
[32 |. A more systematic study of this effect may be a 
promising line of future investigation. 



C. The effect of changing the stress-dependent rate 
function C{s) 

Here we further study the possible effects of details 
of the constitutive behavior on the macroscopic behav- 
ior of the expanding cavity. In this subsection we focus 
on the material function C{s). This phenomenological 
function, as discussed in Sec. IIIBl describes the stress- 



12 



1.5r 



x10"' 



O 



0.5 



0.5 





6000 



FIG. 6: (Color online) The function C{s) of Eq. 1^:23^ with 
C = 7 (dashed line) and of Eq. ^£2^ with A = 30 (solid line). 



dependent STZ transition rates. It is expected to be 
symmetric and to vanish smoothly at s = in ather- 
mal conditions [11| . The plastic rate of deformation for 
s > 1 can be measured in a steady state stress-controlled 
simple shear experiment. For such a configuration the 
deviatoric stress tensor is diagonal and the stable fixed- 
points of Eqs. (|2?20|) - ([2?22l) imply that the steady state 
plastic rate of deformation of Eq. (|2.19p reads 



To 



DP^ = eoe-i/x-C(s) (l-- 



(4.1) 



Therefore, if the steady state relation Xoo{s) is known, 
C{s) can be determined from measuring the steady state 
value of Z)^' for various s > 1 , see for example |2J] . The 
idea then is to interpolate the s — > 0+ behavior to the 
s > 1 behavior with a single parameter that controls 
the amount of sub-yield deformation in the intermediate 
range. In fact, a procedure to measure C{s) at interme- 
diate stresses was proposed in Ref. [l^- Up to now we 
used the one-parameter family of functions ^(s; (^) of Eq. 
(|2.23p . where ( controls the sub-yield deformation. 

We now aim at studying the effect of choosing another 
function C{s). Here we specialize for C{s) = Q{s; A), with 



5(5; A) 



ll + A 



1+ S 



(4.2) 



In Fig. [6] we show C{s) according to the previous choice 
of Eq. (|2.23p with ( = 7 and also C{s) according to the 
present choice of Eq. (|4.2[) with A = 30. The different 
behaviors of C{s) and dC{s)/ds near 5=1 might affect 
differently i?(°-' and R'^^\ thus influencing the stability of 
the expanding cavity. 

In Fig. [7] we compare R^^'' /R^^' as a function of time 
for C{s) of Eq. (|2.23p (previous choice) with C = 7 a-nd 
C{s) of Eq. ft?^ (present choice) with A = 30, both for 
a constant Xoo- An effective temperature perturbation 



FIG. 7: (Color online) i^'^'/i?'"' as a function of time for an 
effective temperature perturbation with x /x =0.03 and 
n = 4 for a°° =4.1. The solid line corresponds to C{s) of Eq. 
(|2.23p with (^ — 7 and the dotted line corresponds to C{s) of 
Eq. (|4.2|l with A = 30, see Fig. [G] In both cases a constant 
Xoo was used. The dotted-dashed line corresponds to C{s) 
of Eq. (|4.2|) with A = 30 and the rate dependent Xoo(toI)^') 
presented in Fig. [T] 



with x'""'^Vx^°'' = 0-03 and n = 4 was introduced at i = 
for 0-°° =4.1. Wc observe that i?'-^'/-^^^^ grows faster for 
the present choice compared to the previous one. For the 
sake of illustration we added the result of a calculation 
with the plastic rate of deformation dependent Xoo as 
discussed in Sec. IIVBI As expected, the effect is mag- 
nified. We conclude that the material-specific function 
of the stress dependence of the STZ transition rates C{s) 
can affect the stability of the expanding cavity, possi- 
bly making it unstable for smaller perturbations. Again, 
the relations between this constitutive property and the 
macroscopic behavior should be further explored in fu- 
ture work. Bringing into consideration explicit macro- 
scopic measurements, one can constrain the various phe- 
nomenological features of the theory of amorphous plas- 
ticity. This philosophy provides a complementary ap- 
proach to obtaining a better microscopic understanding 
of the physical processes involved. 



V. CONCLUDING REMARKS 

We presented in this paper a detailed analysis of the 
linear stability of expanding cavity modes in amorphous 
elasto-viscoplastic solids. The stability analysis is some- 
what delicate due to the non-stationarity of the problem, 
thus a perturbation may grow leaving the problem sta- 
ble if this growth is slower than the growth of the radius 
of the cavity. The radial symmetry of the expanding 
cavity makes it surprisingly resilient to perturbations in 
shape, velocity, external strains and pressure. On the 
other hand the radial symmetry may be lost due to per- 
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turbations in the internal state fields, especially Xj a-nd 
is also sensitive to details of the constitutive relations 
that are employed in the STZ theory. In this respect 
we highlight the role of the plastic rate of deformation 
dependent Xooi'T'nDP'') and of the stress-dependent rates 
of STZ transitions C{s). It is difficult to reach conclusive 
statements, since growth of perturbations beyond the lin- 
ear order invalidate the approach taken here, calling for 
new algorithms involving surface tracking, where the size 
of perturbations is not limited. Nevertheless the results 
point out that instabilities are likely, motivating further 
research into the nonlinear regime. Of particular interest 
is the possibility to select particular forms of constitutive 
relations by comparing the predictions of the theory to 
macroscopic experiments. This appears as a promising 
approach in advancing the STZ theory towards a final 
form. 
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APPENDIX A: DIFFERENTIAL OPERATORS IN 
POLAR COORDINATES 



iv-VT)gg = VrdrTee + —Tre + —Tg,. + —deTe, 

^ '-' /l^ /Y* ff 



(A3) 



Substituting Eqs. (|2.16p for T we obtain the needed ex- 
pressions for [v ■ Ve*^') . . in Eq. (|2.13p . 



APPENDIX B: THE QUASI-STATIC AND 
INCOMPRESSIBLE CASE 

The aim of this appendix is to derive independently the 
linear perturbation theory for a quasi-static and incom- 
pressible case and to compare to the inertial and com- 
pressible case in the limit of large bulk modulus K and 
small velocities v. We show that the results in this limit 
agree, giving us some degree of confidence in the deriva- 
tion and the numerical implementation of the equations 
in both cases. 

The unperturbed problem in the quasi-static and in- 
compressible limit was discussed in detail in [13| and is 
obtained by taking the quasi-static and the incompress- 
ible limits in the equations of Sec. Ill CI Before consider- 
ing the linear stability problem, wc stress that the linear 
perturbation theory of the STZ equations, presented in 
Sec. piIB[) remains unchanged in the present analysis. 
Only the equations of motion and the kinematic equa- 
tions are being modified. In the absence of inertial terms, 
the equations of motion (|2.7p become 



The aim of this Appendix is to derive some expres- 
sions in polar coordinates that were used earlier in the 
paper. Specifically, our goal is to calculate the divergence 
and covariant derivative of a tensor in polar coordinates. 
We represent a general second order tensor T in polar 
coordinates as 



T^T.jBi^ej , 



(Al) 



where e^ and Ej are unit vectors in polar coordinates and 
(8) denotes a tensor product. Using Eqs. (|2.9p - (|2.10|) we 
obtain 



e^ • (V • T) == drTrr 

ee ■ (V • T) == drTre 



Trr , dgTgr T( 



es 



r r r 

Tr9 . Tgr . dgTog 



(A2) 



Substituting Eqs. (|2.1ip for T, we obtain the right-hand- 
sides of Eqs. ((2^2)1 . 

We proceed now to calculate the covariant derivative 
of a tensor v ■ VT. Using Eqs. (j^ . (piU)) and (|XT|) we 
obtain 



ve 



'"e, 



Vg 



Tgr , 



(t;-VT)„ = VrdrTrr + —dgTrr - —T„ 

r r r 

{V ■ VT)^g = VrdrTr-g + —T„ + —OgTrg - —Tgg , 

r r r 



{v ■ VT) 



Br 



1 rr , "e „ , "^e ^ „ vg 

VrOr-Lgr H Jrr H Oglgr Ig, 

r r r 



dr 1 ds 1 dp 2t 
dr r do r do r 



Idr 

r do r^ dr 

To first order we obtain 



I d , y . dp 
r s) 



dr 



inr(i) 



2s(i) 5p(i) as(i) 

= — 1 

r dr dr 



(Bl) 



(B2) 



(B3) 



dr 



+ !!(5(i)_p(i)) + ^ = 0. (B4) 



The boundary conditions of Eqs. (|3.12p - (|3.13|) can be 
further simplified by using the force balance equation to 
zerotli order and the zeroth order boundary conditions 
ofEq. (lOTl) 



dp(o) as(o) 2s 



(0) 



dr dr r 

Substituting into (|3.12p and (|3.13p wc obtain 



(B5) 



,..,«...„,<«,«<..,. a£!!W!l, ,ee, 



'(i?(")) 



2sW(J?W)i?(i) 



(B7) 
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In addition, all the first order fields decay as r^oo. In 
principle, the initial conditions for the partial differential 
equations for the first order fields depend on the type of 
perturbation under consideration. For explicit perturba- 
tions in the shape of the cavity, i.e. i?'^^^(0) ^ 0, we can 
determine the initial stress field by assuming it is sim- 
ply the quasi-static linear elastic solution corresponding 
to the perturbed circle. In order to obtain this solution 
we start with the bi-Laplace equation for the Airy stress 
function $ (sl 



V^V^* = , 
where the stress tensor components are given by 



(B8) 



1 a$ 1 92$ 

r dr ' r^ dO"^ 




92$ d 

dr^ ' ""^^ ~ dr 


/I 9$ 

[r 89 



We then expand the solution in the form 



(B9) 



(BIO) 



The general solutions for $'-'^'(r), that also decay at in- 
finity, are given by 



$(i)(r) 



-ri+2 



+ fer-" , 



(Bll) 



with n>0. Substituting in Eqs. (|B9p . using the bound- 
ary conditions to first order and the following zeroth or- 
der solution 



-nv, 



'Vi- 



(1) 



1 

2^ 

1 
2 

1 

2^ 



r 



dt 



dr 



-+v). 



dr 



^DP 



/(I) 



dvn nVr —Vn 



dr 



(B16) 



(B17) 



dt 



(0)9rW _2£(%^ 
^ dr r 



+ D: 



re 



In order to propagate s^^' and r'^' in time according to 
these equations we need to know Vr and Vg at each 
time step. However, a basic feature of the quasi-static 
problem is that there is no evolution equation for the 

velocity field. Therefore, we must calculate Vr and Vg 
in a different way. 

We now discuss the major mathematical difficulty in 
the quasi-static formulation, i.e. the absence of an ex- 
plicit evolution equation for the velocity field v^^\ To 
overcome this difficulty, we should derive new ordinary 
differential equations for vl and Vg such that their time 
evolution is inherited from the other fields in the prob- 
lem. The first equation can be obtained readily by adding 
(IbTsI) to (IbT6)) 



dv. 



(1) 



vi'^ 



,(1) 



Or 



0, 



(B18) 



from which we can extract v, 



(1) 



.(0) 



•,ee 



one obtains 



It 



(i?(o))' 



'T^^ = , (B12) 



i?(0) 



M^' 



R(0) ' 



i?(o) 



n+2 






(B13) 

The resulting stress components are easily calculated, 
from which we obtain 



,(1) 






s(i) = r(i) = 



RW 



n{l 






n+2' 



(B14) 



These are the initial conditions for the first order stress 
tensor components in terms of the initial R^^'. To pro- 
ceed we expand Eqs. (|2.4p to first order, assuming 
K —i-oo. 



dvi'^ 




dr 
1 
~2jl 


dt 



(B15) 



,,(0) 



9s(i) (i)9s(o) 



dr 



dr 



DP^ 



(1) 




(B19) 



In order to obtain the second equation, we eliminate p^^' 
from the equations by operating with -^^ on Eq. (|B4[) . 
adding the result to Eq. (|B3p and taking the partial time 
derivative to obtain 

9sW 1 8 ( df(^'>\ 2 9f(i) 7if(i) 2s(i) 



dr 



dr 



dr 



dr 



= . 

(B20) 

Here and elsewhere the dot denotes partial time deriva- 
tive. Using Eqs. (|B15|) and (|B17|) we obtain 



^(1) 



^(1) =-' 




(B21) 



= -2^ 



(B22) 



Substituting the last two relations in Eq. (|B20[) and us- 
ing Eq. (|B19p . we obtain a fourth order linear ordinary 
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differential equation for Vr . Since it is straightforward 
to obtain, but very lengthy, we do not write it expUcitly 
here. It is important to note that the coefficients in this 
equation depend on time and therefore by solving it at 
each time step we effectively have a time evolution for the 
velocity field. Once one solves for vi- , Eq. (|B19p can be 
used to calculate Vg . The forth order linear differential 
equation requires four boundary conditions. 

The first boundary condition is obtained by using 
Eq. (Ell, with Eqs. (IbTT]) . (EH) and (|BT91) we ob- 
tain a linear relation between Vr ' (i?^°^), drVr (R'^^^) and 

d'^Vr {R^^'), which is the required boundary condition. 
Another boundary condition is obtained by multiplying 
Eq. (jB4p by r. operating with ^ — dt + vi- dr on the 
result and using Eq. (|B6p . Additional simple manip- 
ulations result in a linear relation between vl (i?'°)), 
drvi^\R^°^), a2„(i)(i?(o)) and d^vi'\R^°^). This is a sec- 
ond boundary relation. Two other boundary conditions 
are obtained from the requirement that Vr vanishes at 
cxD with vanishing derivative 



,(1) 



(oo) = and 



dvr (oo) 



dr 



(B23) 



With these four boundary conditions the forth order dif- 
ferential equation can be solved in the following way: 
at each step we guess Vr\R^'^\t) and drVr '{R^^\t) 
and use the first two boundary conditions to calculate 
d^vi^\R^°\t) and d^vl'^\R^"\t). Then we use the forth 

order differential equation to calculate Vr and drVr at 
CXD. We improve our guess until the solution satisfies Eqs. 
(|B23[) (the shooting method). 

Thus, we have a complete solution procedure (assum- 
ing that the plastic rate of deformation is known, see 
Section IIIIBP : for a given s'^^)(r, i) and T*^^)(r, i) we solve 
the forth order differential equation for Vr ir, t) following 
the procedure described above. Having Vr {r,t) we use 
Eq. ((BT9)) to obtain vl^\r,t). Then we use Eqs. ((BT5| . 
(|BT7)) and ^^ to propagate s^^^r), r'^^^r) and i?(^' in 
time. We follow the same procedure at each time step 
to obtain the full time evolution of the perturbation. We 
note that we have eliminated p^^-* (r, t) from the problem, 
though we can calculate it at every time step using Eq. 
(lB3|) or (lB4l) . 

We are now able to compare the quasi-static and in- 
compressible case to the inertial and compressible coun- 
terpart in the limit of small velocities and large bulk 
modulus K. We introduced at i = a perturbation of 
magnitude /J^^^i?'"' = lO^"'* to the radius of the cavity, 
with a discrete wave- number n = 2, and solved the dy- 
namics in both formulations. We chose a°° < a*'' such 
that the velocities are small and K = 1000 in the iner- 
tial case in order to approach the incompressible limit. 
In Fig. [5] we compare -R^°-' and i?^^-* for both the quasi- 
static and the inertial formulations. The agreement is 
good. Note that the stability of the expanding cavity 
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■- Inertial (K=1 OOP) 




FIG. 8: (Color online) Upper panel: i?'"' as a function of time 
for both the quasi-static (solid line) and the inertial (dashed 
line) formulations. We used a°° = 2, ^ = 15 and K = 1000 in 
the inertial case to ensure small velocities and almost incom- 
pressibility. Lower panel: W^' as a function of time for both 
the quasi-static (solid line) and the inertial (dashed line) for- 
mulations with a discrete wave-number n = 2. The agreement 
of the curves in both panels is favorable, where the small dif- 
ferences are attributed to the finite bulk modulus K in the 
inertial case. 
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FIG. 9: (Color online) The fields s''°'>{r), s'-^\r), x^°\r) and 
X (t) at i = 2500to for both the quasi-static (solid line) and 
the inertial (dashed line) formulations. The parameters are 
the same as those of Fig. [S] 



depends on the time dependence of the ratio _R'^^'/i?'"^; 
however, we do not discuss the stability yet, but focus on 
the comparison between the two formulations. In Fig. [9] 
we further compare the predictions of the two formula- 
tions for the zeroth and first order deviatoric stress field 
s and effective disorder temperature x at a given time. 
In all cases the differences are practically indistinguish- 
able. We thus conclude that the quasi-static formulation 
and the inertial one agree with one another, giving us 
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some confidence in the validity of both. In particular we 
conclude that the inertial formulation can be used with 



confidence also for high velocities where the quasi-static 
counterpart becomes invalid. 
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